Neutron matter based on consistently evolved chiral three-nucleon interactions 
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We present the first results for the neutron matter equation of state (EOS) using nucleon-nucleon 
and three-nucleon chiral effective field theory interactions that are consistently evolved in the frame- 
work of the Similarity Renormalization Group (SRG). The dependence of the EOS on the SRG 
resolution scale is greatly reduced when induced three-nucleon forces (3NF) are included and the 
residual variation, which in part is from missing induced four-body interactions, is comparable to 
estimated many-body perturbation theory truncation errors. The relative growth with decreasing 
resolution of the 3NF contributions to the energy per neutron is of natural size, but it accelerates at 
the lowest resolutions where strong renormalization of the long-range 3NF matrix elements is also 
observed. 
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Chiral effective field theory (EFT) [1 offers a sys- 
tematic expansion of nuclear forces well suited to meet 
the calculational challenges of neutron matter, which 
span the extremes of low-density universal properties 
to the dense matter in neutron stars. Conversely, neu- 
tron matter provides a powerful laboratory for testing 
chiral EFT power counting at relevant nuclear densi- 
ties, since only long-range three-nucleon forces (3NF) 
contribute at next-to-next-to- leading order (N^LO) [2j 
and there are no new parameters for three-nucleon (3N) 
and four-nucleon interactions at next-to-next-to-next-to- 
leading order (N^LO) [1 . For example, it was recently 
shown that 3NF at N^LO give relatively large contribu- 
tions to the neutron matter equation of state (EOS) [3l|4], 
which may indicate that a chiral EFT with explicit delta 
degrees of freedom would be more efficient. 

At present, the largest uncertainties in microscopic cal- 
culations of neutron matter based on chiral EFT interac- 
tions are because the low energy constants in the Hamil- 
tonian are not all well determined. This leads to uncer- 
tainties in observables such as the nuclear symmetry en- 
ergy and radii of neutron stars [5 . However, direct calcu- 
lations based on chiral interactions using many-body per- 
turbation theory (MBPT) also have non-negligible the- 
oretical uncertainties due to truncations of the many- 
body expansion [3l|4]. Renormalization group (RG) evo- 
lution of nuclear interactions to lower resolution scales 
significantly improves the convergence of MBPT, but 
in prior calculations three-nucleon interactions have not 
been evolved consistently [2 . Here we present the first 
results for the neutron matter equation of state based on 
consistently evolved chiral nucleon-nucleon (NN) and 3N 
forces (see Fig. [T]). These results show how RG transfor- 
mations can enable simplified and efficient many-body 
calculations for neutron matter with controlled theoreti- 
cal error bars. 

We build upon Ref. |6^, which presented a framework 
for the simultaneous evolution of nuclear NN and 3N in- 
teractions in a continuous plane-wave (momentum) basis 
via the Similarity Renormalization Group (SRG). This is 



an alternative to using a harmonic oscillator basis [7^ . It 
provides independent checks of the SRG evolution, easier 
access to alternative SRG generators (see Refs. [S'.'Q^) and 
a means to test approximations for induced 3NF [6 . The 
momentum-space matrix elements can be easily trans- 
formed to an oscillator basis for use in calculations of 
finite nuclei by configuration interaction [10^, coupled 
cluster [mill], in-medium SRG methods [13 [12, self- 
consistent Gorkov Green's function theory [T5^ or for nu- 
clear shell-model calculations [16 . 

Most important in the present context is that these 
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FIG. 1. (Color online) Energy per neutron as a function of 
neutron density for different SRG resolution scales. The re- 
sults are grouped according to whether no induced 3NF are 
included (NN-only), the induced 3NF are included but no 
initial 3NF (3N-induced), or initial and induced 3NF are in- 
cluded (3N-fun). The initial interaction is the 500 MeV N^LO 
NN potential from Ref. [IT] combined with the N^LO 3NF 
using the consistent low-energy constants ci = — 0.81GeV~^ 
and C3 = -3.2 GeV"^ [2]. 
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FIG. 2. (Color online) Partial- wave convergence of the Hartree-Fock 3NF-only energy per neutron as a function of density at 
two different resolution scales. The left panel shows the results before evolution (A = oo) and the right panel at A = 2.0 fm~^. 
The Hamiltonian is the same as for the 3N-full results in Fig.[l] For A = oo, only J < 5/2 give significant contributions so the 
lines for the higher partial waves are nearly indistinguishable from the exact result. 



evolved interactions can be directly applied in micro- 
scopic calculations of the nuclear EOS. Applications to 
the triton in Ref. [6 verified that nuclear interactions 
at low resolution scales contain much weaker couplings 
between low and high-momentum states. This renders 
the nuclear many-body problem more perturbative and 
therefore more tractable, with all low-energy observables 
preserved exactly when all induced contributions are in- 
cluded. However, recent results for finite nuclei using 
oscillator-evolved 3NF found that high-momentum parts 
of the chiral two-pion-exchange 3NF led to significant 
resolution-scale dependence and overbinding in medium- 
mass nuclei [T8^, HlQ]. Neutron matter provides a test 
laboratory (although limited, because only neutrons) to 
study if these systematics carry over to infinite matter in 
MBPT with the new evolution scheme. We have direct 
access to the form of the 3NF and the scaling behavior 
of their contributions to the energy as a function of the 
resolution scale. 

The SRG flow equations we solve can be written in the 
form [20] 

where Hg = T^ei + Vg denotes the Hamiltonian as a func- 
tion of the flow scale parameter 5, and r]s labels the gener- 
ator of the KG transformations. In practice it is more in- 
formative to replace s with the resolution (or decoupling) 
scale A = which has units of momentum. Here we 

choose r]s = [^rei, Hs] with the relative kinetic energy Trei, 
as in Ref. [6] and most prior investigations. With this r]s 
the flow equation generates a continuous series of unitary 
transformations that renormalizes the Hamiltonian (and 
all other operators), driving Hg towards a diagonal form 



in momentum space [2j . We recast Eq. ([T]) into separate 
flow equations for the matrix elements of the NN and 3N 
interactions [6l [20] and solve them simultaneously in a 
momentum partial- wave basis. 

For the NN forces we use a standard partial- wave basis 
of the form |p; {LS)JT)^ where p is the relative momen- 
tum and L, 5, J and T are the orbital angular momen- 
tum, spin, total angular momentum and isospin of the 
interacting pair. For the three-body basis we choose 

Iftft; [{LS)J{ls,)j] JMTU)T%) = \pqa)^ , (2) 

where pi and Qi are the three-body Jacobi momenta of 
particle i. The quantum numbers /, Si = 1/2, j and 
ti = 1/2 are the orbital angular momentum, spin, to- 
tal angular momentum and isospin of particle i relative 
to the center-of-mass of the pair with momentum p. J 
and T are the total three-body angular momentum and 
isospin quantum numbers (for details see Refs. [22} |23]). 
The 3NF at chiral order N^LO are independent of the 
projections and Tz [T- We use a to abbreviate the 
angular momentum and isospin quantum numbers. 

In the SRG evolution we take all NN interaction ma- 
trix elements up through Jmax = 7 into account, and for 
the 3NF matrix elements we include all interactions up 
through J 9/2 and J^ax ^max 5. In Fig. [5] we 
show the convergence of the partial-wave 3NF contribu- 
tions in the Hartree-Fock approximation as a function of 
density for the initial interaction and after evolution to 
A = 2.0 fm~^. While the convergence pattern is some- 
what altered after evolution, the partial-wave truncation 
is reliable at both scales and the results are well con- 
verged. The exact Hartree-Fock energy can be calcu- 
lated without a partial wave expansion for the unevolved 
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interaction [2l [24] , and in this case the energy through 
J' = 9/2 is converged within 0.4 percent at nuclear sat- 
uration density. 

The RG evolution of the 3NF can be performed in- 
dependently for different values of J7, T and three-body 
parity tts = (— 1)^+^ For applications to neutron mat- 
ter we only need matrix elements for the isospin channel 
T = 3/2. Our calculations are based on antisymmetrized 
matrix elements of the form 



where v/23 is the i-th Faddeev component of the three- 
body interaction, ^123 = (1 + ^123 + ^132) is an antisym- 
metrizer, and ^123(^132) is the cyclic (anti-cyclic) per- 
mutation operator (see Ref. ^2). Such antisymmetrized 
matrix elements with two antisymmetrizers, in contrast 
to partially antisymmetrized matrix elements typically 
used in Faddeev calculations (see, e.g., Ref [25]), are 
particularly suitable for solving the SRG flow equations 
for 3NF (see Ref. [6]). The initial matrix elements are 
generated directly in this form using a novel automa- 
tized partial- wave decomposition for the 3NF [25l [26] . 
This new method ensures that the interaction is anti- 
symmetrized exactly in spin and isospin space. 

Due to the softening of NN and 3N forces during the 
RG evolution, it is possible to apply MBPT for the calcu- 
lation of the neutron matter equation of state at low res- 
olution scales. In the present calculations we resum con- 
tributions from NN forces in the ladder approximation 
while contributions from 3NF are calculated in Hartree- 
Fock approximation. Based on results of Ref. [2 we ex- 
pect this approximation to provide the dominant 3NF 
contributions at low-resolution scales. We also compute 
the NN second-order contributions for comparison, which 
allows us to probe the perturbativeness of the NN inter- 
actions as a function of the RG scale A. The inclusion 
of higher-order diagrams involving 3NF requires signif- 
icant computational storage because of the coupling of 
partial waves with different J and tt. Such calculations 
are currently in progress. 

In the Hartree-Fock approximation, the 3NF contribu- 
tions to the energy per volume is given by: 



momentum basis can be expressed in the form 

^ (pq^M^i^|^123'1^123A23|pq5'M^^) 
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where nk are the zero-temperature occupation numbers. 
Compared to the relations of Refs. [21 [2^ an extra fac- 
tor 1/3 appears because the 3NF are antisymmetrized in 
the initial and final states. The spin sum in the Jacobi 
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= 7iE^««' E Stj^JjJ'j'LL'W 
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using standard notation for the angular momentum cou- 
pling, the Clebsch-Gordan coefficients CjjJ^^^^^^, and 

X yj2x + 1. 

In Fig.[l] we show the results for the energy per neutron 
at four different resolution scales A as a function of neu- 
tron number density n. This range of resolution scales 
has also been used in in previous studies of neutron mat- 
ter and nuclear matter [21^ based on l^ow /c- evolved NN 
interactions. The present calculations use the N^LO NN 
potential (A = 500 MeV) of Ref. [17] plus the consistent 
3NF at N^LO with the couplings ci = -0.81 GeV"^ and 
C3 = — 3.2GeV~^ at the initial resolution scale A = 00 
(s = 0). The other components of the fuh N^LO 3NF give 
no contributions in neutron matter [2]. We calculate the 
EOS in three ways: NN-only, 3N-induced and 3N-full (see 
Ref. [7 and the caption to Fig.jl]). The NN contributions 
are resummed and the 3NF contributions are calculated 
in Hartree-Fock approximation. When induced 3NF are 
taken into account, we find a dramatically reduced A de- 
pendence over the entire density range compared to in- 
cluding only NN forces. When initial 3NF are also in- 
cluded the spread of the results increases, but remains 
still significantly smaller than the spread of the NN-only 
results. The 3N-full energy per neutron for A = 2.0 fm~^ 
at saturation density is about 0.5-1.5 MeV higher than 
found in calculations with the same initial NN interac- 
tion but with the unevolved N^LO 3NF included with a 
cutoff ranging from 2.0-2.5 fm~^ 0. More detailed com- 
parisons will be made in a future publication. 

In Fig. [3] we show the energy per neutron at nuclear 
saturation density for two different NN interactions as 
a function of the resolution scale. The left panel uses 
the NN and 3N interactions as in Fig. [l] and the right 
panel uses the N^LO potential (A/A = 450/500 MeV) 
of Ref. [29] plus the consistent 3NF with the low-energy 
constants ci = -0.81 GeV"^ and C3 = -3.4 GeV"^ The 
bands in both panels at the right side indicate the size of 
unevolved second-order contributions containing NN and 
3N forces; these diagrams can be evaluated at the initial 
RG scale using the framework of Ref. 0. The band is 
smaller in the right panel due to the smaller initial cutoff 
of the NN potential. 

The variation of the energy in the 3N-full case is al- 
ways within the width of the band. This suggests that the 
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FIG. 3. (Color online) Energy per neutron as a function of the resolution scale A at nuclear matter saturation density 
Us = 0.16 fm~^ for two different interactions. The left panel uses the NN and 3N potentials as in Fig. [l] and the right panel uses 
the N^LO potential (A/A = 450/500 MeV) of Ref. [29] plus the consistent 3NF with the low-energy constants a = -0.81 GeV~^ 
and C3 = — 3.4GeV~^. The same three calculations (NN-only, 3N-induced, and 3N-full) from Fig. [l] are used here. The solid 
lines are energies using a resummed particle ladder sum of the NN contributions, with the straight dotted lines marking the 
energies at A = 3.0 fm~^. The dashed lines are energies including up to second-order for the NN contributions. The shaded 
areas at the right indicate the magnitude of the attractive second-order contribution from unevolved NN and 3N forces. 



inclusion at finite A of neglected second-order diagrams 
with 3NF, which are attractive and should decrease in 
magnitude with decreasing A [2], may systematically re- 
duce the observed variation with A above 2.0 fm~^. The 
small increase in energy below A = 2.0 fm"^ for both 
3N-full and 3N-induced curves might be attributed to 
induced four-body forces, but it would be premature to 
draw quantitative conclusions. 

The inclusion of induced 3NF contributions greatly re- 
duces the resolution-scale dependence of the EOS even 
at the present truncation of the many-body expansion. 
In the left/right panel of Fig. [s] we find a maximal energy 
variation of about 390/470 keV for the 3N-induced cal- 
culations and 650/450 keV for the 3N-full calculations at 
saturation density. In comparison, we find a total vari- 
ation of about 3.6/4.2 MeV when 3NF are completely 
neglected. Furthermore, Fig. [3] also demonstrates the 
increased perturbativeness of the many-body expansion. 
The solid lines show results using the NN ladder sum 
while the dashed lines are use diagrams up to second or- 
der. For the NN interaction at N^LO (right panel), dia- 
grams beyond second order in the particle-particle chan- 
nel give only very small contributions. 

The net relative growth of two- and three-body con- 
tributions to the energy at two densities are shown in 
Fig. [4j with and without initial 3NF. The solid lines are 
for the N^LO 500 MeV NN potential while the dashed 
lines are for the N^LO 450/500 MeV potential. The size 



of the initial 3NF sets the scale of a natural ratio at that 
density. It is evident that the change of the net ratio 
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FIG. 4. (Color online) Scaling of the contributions from NN 
and 3N forces as a function of A. The dashed lines show the 
ratio of these two contributions based on the 500 MeV N^LO 
NN potential 17 and the dashed lines based on the 450/500 
MeV N^LO potential ^ (see also Fig.jsJ. 
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FIG. 5. (Color online) Matrix elements of the evolved 3N potential (^a = l|V^i23|^^<^^ — l) (see main text) for J = 1/2, 
T = 3/2 and positive total parity tts at the hyperangle = tt/S (bottom). The interactions are the same as in Fig. IT] 



with A remains natural even down to the smallest values. 
There is no obvious trend with density; such a trend may 
be obscured by cancellations among contributions to the 
net energies. 

We can also examine the evolution of the 3NF by look- 
ing at slices of the matrix elements. To do this we intro- 
duce the hyperradius = -j-S/Aq'^ and the hyperangle 
tan^ = 2p/{\^q) and visualize the matrix elements as a 
function of ^ at a fixed hyperangle. A representative ex- 
ample is shown in Fig. [s] for = tt/S for the dominant 
partial wave with a = and L = J = 6' = / = 0, T = l, 
j = 1/2, J = 1/2 and tts = 1. In this case we observe 
softening from A = ootoA = 2.0 fm~^. Further evolu- 
tion causes a strongly attractive part to appear at small 
momenta; however, its impact will be mitigated by phase 
space factors. This renormalization merits further study. 

In this paper we have presented the first neutron mat- 
ter calculations based on a fully consistent RG evolution 
of two- and three-body chiral EFT interactions. Includ- 
ing induced 3NF greatly reduces the resolution scale de- 
pendence of the neutron EOS compared to NN-only cal- 
culations, with the residual dependence comparable to 
the expected magnitude of omitted many-body correc- 
tions. Thus we are not able to make a definitive state- 
ment about the size of induced four-body contributions, 
but there are no indications of unnatural growth. Future 
calculations will include both neutron and nuclear mat- 
ter with higher-order (beyond Hartree-Fock) diagrams 
including 3NF, which will allow a more complete assess- 
ment of higher-body contributions. We will also calculate 
error bands based on uncertainties in the 3NF input. In 
the future it will be also straightforward to include contri- 
butions from 3NF at N'^LO [30l|3T] once initial matrix el- 
ements in partial wave representation are available. This 
work is currently in progress. 
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